Cite this Article
Xuan Li, He Bin, Hu Li-Fa, Li Da-Yu, Xu Huan-Yu, Zhang Xing-Yun, Wang Shao-Xin, Wang Yu-Kun, Yang Cheng-Liang, Cao Zhao-Liang, Mu Quan-Quan, Lu Xing-Hai. Configuration optimization of laser guide stars and wavefront correctors for multi-conjugation adaptive optics. Chinese Physics B, 2016, 25(9): 094216  
Permissions
Configuration optimization of laser guide stars and wavefront correctors for multi-conjugation adaptive optics
Xuan Li1, He Bin1, 2, Hu Li-Fa1, †, , Li Da-Yu1, Xu Huan-Yu1, Zhang Xing-Yun1, Wang Shao-Xin1, Wang Yu-Kun1, Yang Cheng-Liang1, Cao Zhao-Liang1, Mu Quan-Quan1, Lu Xing-Hai1
State Key Laboratory of Applied Optics, Changchun Institute of Optics, Fine Mechanics and Physics, Chinese Academy of Sciences, Changchun 130033, China
University of Chinese Academy of Sciences, Beijing 100049, China

 

† Corresponding author. E-mail: hulifa@ciomp.ac.cn

Project supported by the National Natural Science Foundation of China (Grant Nos. 11174274, 11174279, 61205021, 11204299, 61475152, and 61405194) and the State Key Laboratory of Applied Optics, Changchun Institute of Optics, Fine Mechanics and Physics, Chinese Academy of Sciences.

Abstract
Abstract

Multi-conjugation adaptive optics (MCAOs) have been investigated and used in the large aperture optical telescopes for high-resolution imaging with large field of view (FOV). The atmospheric tomographic phase reconstruction and projection of three-dimensional turbulence volume onto wavefront correctors, such as deformable mirrors (DMs) or liquid crystal wavefront correctors (LCWCs), is a very important step in the data processing of an MCAO’s controller. In this paper, a method according to the wavefront reconstruction performance of MCAO is presented to evaluate the optimized configuration of multi laser guide stars (LGSs) and the reasonable conjugation heights of LCWCs. Analytical formulations are derived for the different configurations and are used to generate optimized parameters for MCAO. Several examples are given to demonstrate our LGSs configuration optimization method. Compared with traditional methods, our method has minimum wavefront tomographic error, which will be helpful to get higher imaging resolution at large FOV in MCAO.

PACS: 42.70.Df;42.30.Wb;42.68.Bz;42.15.Fr
Keyword:laser guide star;liquid crystal wavefront corrector;adaptive optics;phase reconstruction
1. Introduction

Because of atmospheric random disturbances, the wavefront from a star will be distorted, which blurs the image captured by large aperture, ground-based optical telescopes. Adaptive optics (AOs)[1,2] technology can be used to restore the imaging resolution close to the telescope diffraction limit by compensating the aberration induced by atmospheric turbulence in real time. The AO system is now considered to be an essential setup of many new large aperture optical telescopes.[3,4] In traditional adaptive optics, the distortion is measured by a wavefront sensor (WFS) and corrected by a wavefront corrector (WFC), such as a deformable mirror (DM)[5] or liquid crystal wavefront corrector (LCWC),[6–13] conjugated to the telescope pupil.[14,15] Unfortunately, this approach only achieves a good imaging resolution within a limited field of view (FOV) smaller than the isoplanatic angle that is generally only a few arc-seconds at visible wavelengths. Anisoplanatism leads to a limited FOV, which is a bottleneck of many applications of AO systems. To enlarge the FOV, multi-conjugate AO[16,17] has been developed by using several natural guide stars (NGS) or laser guide stars (LGS) in different directions. In a multi-conjugation adaptive optics (MCAOs) system, such as GeMS,[18,19] every guide star is sensed by an independent WFS, and several WFCs are conjugated to turbulence layers at different heights.

The atmospheric tomographic phase reconstruction and projecting three-dimensional turbulence phase onto WFCs are a necessary and very important calculation step for MCAO, which has been investigated by many groups.[20–24] These groups have discussed the influence of the number of mirrors and LGSs, but did not consider optimizing the specific location of LGSs and WFCs further. In these reported works, LGSs were simply located at the vertices of inscribed regular polygon inside the FOV, and the WFCs were placed at a position conjugated to the weighted mean height. Although Fusco proposed a method to optimize the NGS separation for a given FOV in Ref. [25], the distribution of NGS is not controllable, unlike the LGS that can be controlled by laser beam launch telescope. Up to now, LGS distribution and reasonable WFC conjugation altitude have not been optimized, though they are very important to obtain a high accuracy of atmospheric tomography. To resolve this problem and improve system imaging resolution at large FOV for MCAO, we have presented a new method to optimize the LGSs’ configuration and LCWCs’ conjugation heights. Our method is based on evaluating the field-averaged mean-square residual phase error variance of MCAO across the whole FOV for arbitrary configuration of LGS and LCWC analytically based on Zernike modes. The error variance is decomposed into reconstruction error related to LGS and correction error related to LCWC. Then, on one hand, by minimizing the wavefront reconstruction error averaged over the FOV, we get the optimal LGS configuration. On the other hand, we get the optimal conjugate altitudes of LCWCs by minimizing the error of wavefront correction for the MCAO system. After the optimization of the location of LGSs and LCWCs, we can obtain a phase reconstruction and compensation with high precisions, which is of great significance for multi-LGS AO system design.

In Section 2.1, Zernike coefficients representation for MCAO are given. Using a modal decomposition of the phase onto the Zernike polynomials basis, a simple analytical solution is proposed. The equations of shape of LCWCs and correspondent performance of MCAO are derived in Section 2.2. In Section 2.3, we also derive formulas for the computation of the covariance matrices based on the statistical property of atmosphere turbulence. In Section 3, a demonstration of multi-LGS optimization process is given and discussed. Finally, our conclusions are given in Section 4.

2. Theory
2.1. Zernike coefficients representation for MCAO

Modal decomposition of wavefront phase using Zernike polynomials is widely used in theoretical studies of AO. For a unit circular aperture without obstruction, the Noll Zernike modes are defined as follows:[26]

where

The index n and m are called the radial degree and the azimuthal frequency, respectively. Any distorted wavefront can be decomposed based on Zernike modes on its unit circular. Note that each LCWC is assumed to correct spatial frequencies up to a cut-off frequency and the reconstruction modes of each WFS is limited by the number of sub-apertures. We assume that the series representation of all phase error in this paper is truncated with finite maximum order J. The minimum order is 4 without considering piston and the global tip-tilt component that is detected by NGS and compensated by tilt mirror.[27,28]

Figure 1(a) shows the simple scheme that light beam from a star at direction α within the telescope’s FOV propagates through a turbulence layer conjugated to LCWC at height h and enters the telescope’s pupil. As shown in Fig. 1(a), the meta-pupil is the large circle with a diameter Dh = Dp + hθfov at height h. Atmospheric turbulence aberration flowing through the larger circle as shown in Fig. 1(b) at height h will affect the imaging resolution of a sky object seen by the telescope. The Zernike basis on this circular pupil is denoted as {Z(2r/Dh)}. The light beam from the target illuminates only a portion of the meta-pupil, as shown by the gray area in Fig. 1(b). The beam footprint is enclosed by a small circle with a diameter Dp. In addition, the displacement between the center of beam footprints and that of the meta-pupil is equal to hα. The basis on the small-circular pupil is {Z(2r/Dp)}. Tokovinin[29] gave a mode projection matrix P which permits to represent a portion of a large-circle mode by a sum of small-circle modes

where Pji stands for the relationship between the i-th Zernike mode in meta-pupil and the j-th Zernike mode in small pupil. This can be calculated by

For certain Dp and θfov, it is dependent on α, h. The last written will be simpler using a Cartesian polynomial representation.[30]

Fig. 1. Scheme of mode projection, Dp is the diameter of telescope, θfov is the FOV.

It is reasonable to assume that atmospheric turbulence aberration can be optically separated into L layers in different heights hl (l = 1,2,…,L)

where al is the Zernike coefficients for the l-th turbulence layer. It is a vector with (J-3) elements

where al, j is the coefficients associated of the j-th Zernike mode of the l-th turbulence layer. The radius ρl of a meta-pupil for the l-th turbulence is defined as

Let us consider a target star at direction αt. Light from star in direction αt will pass through a series of atmosphere layers in small circles at different heights into the telescope. In the near-field approximation, the aberration to be corrected on the telescope pupil is

The piston mode is not included because it does not affect the imaging resolution. The coefficient ci is

If there are M LCWCs conjugated to different altitudes h̃m, m = 1, …, M, the surface shape of the LCWCs to compensate the turbulence could be described by

where ãm is a vector of Zernike coefficients for the surface shape of the m-th LCWC which is defined on a meta-pupil with the radius ρ̃m = Dp/2 + h̃mθfov/2. It is a vector with (J-3) elements that piston and tip-tilt are removed

Therefore, the turbulence aberration at direction αt to be compensated by LCWC could be expressed as

where the coefficient c̃i is

According to the orthonormality property of the Zernike polynomials on the unit circle, the instantaneous value of the field-averaged mean-square residual phase error variance could be given by

where x̄ is the average value over the FOV. PTL is the projection matrix of turbulence

and PLC is the projection matrix of LCWC

where

It is obvious that residual wavefront error after correction also depends on the location and shape of LCWCs.

2.2. Compute the shape of LCWCs and performance estimation

This section will discuss how to compute the instantaneous shape of LCWCs that are specified as the coefficients of wavefront expansion on Zernike modes. In Fig. 2, a scheme of multiple-LGS system is given. Light beams from several LGSs are used to measure atmospheric tomographic aberration, and only one LGS is illustrated in the figure. Assume that we have G different LGSs located at height hg in directions βg(θg, ϕg). θg is the zenith angle, and ϕg is the azimuth angle (g = 1,2, …, G). The light beam from LGS origins from a point at finite height and propagates through a cone as shown in Fig. 2. In the near-field approximation, to a position of r on the telescope pupil, the turbulence-induced wavefront phase incoming wavefront phase Λg(r) associated with each LGS is given by

where Φ(l) is the phase funciton of the l-th turblence layer and bg, j is series expansion coefficients of Λg(r). is the Zernike coefficients of the turbulence aberrated wavefront in the cross section of the l-th layer and given by

Fig. 2. Scheme of turbulence tomography using LGS. The gray disks are the footprints of LGS light cone at the correspondent turbulence layer, respectively. Their radius is proportional to a factor Kl = (hg − hl)/hg. The position vector of the center at the l layer is βghl. Light beam from LGS g goes to a position of r through the point Klr + βghl at the l layer.

The number of wavefront reconstruction Zernike modes for WFS is limited by its number of sub-apertures. It is assumed that the highest order of WFS is equal to that of LCWC. Hence, the detection signal including all LGSs is

where b is a column vector with G sub-vectors, each having (J − 3) elements.

We neglect the temporal aspects of the MCAO operation by supposing that all measurements and corrections are done instantaneously. The computation of the optimal mirror shapes involves two sub-problems: atmospheric reconstruction and the derivation of the optimal mirror shape. We converge the two steps at a control matrix R, and the LCWC actuator command vector ã is calculated by

The optimal modal wave-front compensation has a minimization of the statistical expectation of . Therefore, we have to minimize the error to find the matrix R. The expected residual error is quadratic in matrix R. Thus, the value of R that minimizes subject can be determined by partial differential, according to Eqs. (14) and (21)

where the operator 〈•〉 represents the statistical averaging. From Eq. (22), an optimal control matrix R is given by

where Cxy is the covariance matrices of column vector x with column vector y

Putting Eqs. (21) and (23) into Eq. (14) leads to

where Tr stands for the matrix trace operation that is the sum of the diagonal elements. For a certain turbulence model, Dp and θfov, depend on PLC, Cab and Cbb, that are all related to configuration of LGS and LCWC. In order to optimize configuration conveniently, we divide into three parts as follows:

The error of reconstruction JGS is only related to Cab and Cbb that depend on the configuration of LGSs. The error of correction JLC is only related to PLC that depends on the configuration of LCWCs. ΔJ is the common part of JGS and JLC, which is much smaller than JGS and JLC. Thus, we can get the optimal LGS configuration and conjugate altitudes of LCWCs by minimizing JGS and JLC, respectively.

2.3. Computation of the covariance matrices

The computation of the covariance matrices Caa, Cab, and Cbb is necessary for performance estimation and configuration optimization. In this section, we derive formulas for the computation of the three matrices based on the statistical property of atmosphere turbulence.

According to Eq. (19), the turbulence layers are taken to be statistically independent,

These covariances can be computed by a priori information on turbulence and the geometry of the LGSs of AO systems. Whiteley et al. introduced a generalized analysis geometry and used this aperture-and-source geometry with conventional methods to get a general expression for the inter-aperture cross correlation of the Zernike coefficients.[31] Assuming each turbulence layer’s statistics follows a von Kármán power spectral density with outer scale L0 and ignoring inner scale, we can get the following formula:

where BA1i,A2j,l in the turbulence layer l is the covariance between Zernike mode i measured within aperture A1 from the source S1 and Zernike mode j measured within aperture A2 from the source S2. R1l and R2l are the radius of the source light cones at the intersection with turbulent layer l to the corresponding pupil radii at the aperture. pl is the displacement vector between the center of the two intersection. ni, mi, and nj, mj are the radial and the azimuthal orders of the i-th and the j-th Zernike modes, respectively. Jλ(•) is the Bessel function of the first kind of order λ. δ is the Kronecker delta. x0 = πDp/L0. r0 is the Fried parameter of the turbulence. wl is the relative strength of a layer in the atmosphere with

This is dependent on the refractive-index structure constant profile.

The aperture and source geometry for single-layer turbulence is shown in Fig. 3. We can get correspondent R1l, R2l, and pl, according to their A1, A2, S1 and S2, respectively. Now we have three kinds of covariance as

If the Zernike model number is enough, we can also compute covariance matrices approximately by the linear method.[23]

Fig. 3. Analysis geometry used for computing the cross correlation of Zernike coefficients. al is the same as a vector of Zernike coeffcients associated to an infinite altitude beacon in the center of the FOV, measured within a virtual aperture with radius ρl centered at the origin which is indicated with dashed circle. The solid circle at the bottom is the telescope aperture. (a) Covariance between ai,j with ; (b) covariance between with ; (c) covariance between al,i with al,j.
3. Results and discussion
3.1. Turbulence model and system parameters

For this work, the well-known Hufnagel–Valley 5/7 model is used to describe the profile. Troxel[32] mentioned that four layers could be used to describe accurately the turbulence. The altitudes and strength of the turbulence layers are listed in Table 1. The Fried parameter for the whole atmosphere is r0 = 0.12 m at λ = 500 nm, and the outer scale is infinite.

Table 1.

Four-layer turbulence model altitudes and relative layer strength.

.

We consider an 8-m telescope with several sodium LGSs (hg = 90 km). The LGSs are positioned in a circle of radius θg. In this paper, we focus our work on the optimization of configuration of LGS and conjugated heights of LCWCs. For simplification, the fitting error (variance of modes higer than J) and noise modeling are not considered here. Hence, we consider low-order correcting Zernike mode numbers 4-136 (radial order 15) without noise. For our atmospheric model, the uncorrected variance of those modes is 142.5 rad2.

3.2. Optimization of LGS configuration

We compute the reconstruction error of different LGS constellations. The number of LGSs is in the range from 2 to 5. All LGSs are located at the vertices of regular polygons as shown in Fig. 4. On the basis of symmetrical characteristic, different azimuth angles of LGSs are equivalent. In order to optimize the LGS configuration we have to find an optimized θg.

Fig. 4. Illustration of the LGS constellation. (a) G = 2, (b) G = 3, (c) G = 4, (d) G = 5.

The result is shown in Fig. 5, where the error decreases with the increase in the number of LGSs. However, further increase in LGSs’ number will lead to limited decrease of residual error. Therefore, there is a tradeoff between residual error and AOS cost caused by increasing the number of LGSs. For every kind of LGS arrangement, the reconstruction error JGS has a minimum value at a specific θg, which results in an optimal θg. The larger θg is, the weaker the correlation between turbulence layers and the detected signal will be. The smaller θg is, the stronger the correlation of measurement by different LGSs will be, which leads the efficiency of the measurement to reduce because of over-sampling on the turbulence. This indicates that the optimal θg is not the same as usually used. The wavefront reconstruction errors JGS at optimal θg and half FOV 30″ are listed in Table 2, respectively.

Fig. 5. JGS as a function of θg. Gray dashed line at 30 arcsec is the inscribed θg. Results shown are for hg = 90 km, θfov = 1 arcmin, Dp = 8 m.
Table 2.

The optimal θg and corresponding minimum JGS.

.

As shown in Fig. 6(a), for five LGSs and a telescope with diameter of 8 m, the optimal θg increases monotonously as FOV becomes larger. We calculate a series of optimal θg for θfov in the range from 0 to 120 arcsec. For comparison, five LGSs are arranged uniformly at the edge of FOV according to the traditional methods,[21,22,24] and corresponding results are shown by the square line in Fig. 6(a). The θg of both methods varies nearly linearly with θfov, but has different slopes. They are the same for 42″ FOV. Their corresponding minimum average variances are also calculated. The wavefront reconstruction errors for the two sceneries are shown in Fig. 6(b). This indicates that our optimal arrangement of LGSs is obviously better than traditional one, and its reconstruction error is much less than the traditional method for FOV smaller or larger than 42″.

Fig. 6. Comparison of simulation results of θg and JGS with our optimal result (circle) and traditional method that inscribed in the FOV (square) for hg = 90 km, G = 5, θfov = 1 arcmin, Dp = 8 m: (a) θg as a function of θfov; (b) JGS as a function of θfov. Square line: traditional method; circle line: optimal method.
3.3. Optimization of LCWC configuration

Considering the cost and complexity of MCAOS, one generally uses several LCWCs in MCAOS. Here, in our simulation, we use three LCWCs in our simulation with one LCWC conjugated to the ground turbulence layer at 0 km as shown in Fig. 7(a). LCWC1, LCWC2, and LCWC3 are conjugated to turbulence layer 1, 2, and 3, respectively. The simulated results are shown in Fig. 7(b). Especially, LCWC1 coincides with LCWC2 for h̃2 = 0 km, which is equivalent to the case of two LCWCs with one LCWC conjugated to 0 km. This indicates that the optimal height is [0 km, 10 km] for two LCWCs and [0 km, 2 km, 11 km] for three LCWCs, respectively.

Fig. 7. Optimization of LCWCs’ conjugation heights. (a) Scheme of three LCWCs and their conjugated heights. (b) The optimal result of LCWC for 1-arcmin FOV.

Figure 8(a) shows that the reconstruction error JLC for two LCWCs with one at 0 km and another at different heights varying with different FOVs. The curves from the bottom up vary in FOV from 0.2 arcmin to 2 arcmin with the increment of 0.2 arcmin. This indicates that the optimal height conjugated to the second LCWC is always 10 km for different FOVs. The result is not affected by the FOV. The simulation results for three LCWCs are also obtained, but not shown here. Similarly, a series of optimized conjugation heights at different FOVs are also obtained. As a comparison, figure 8(b) shows that the minimum reconstruction error JLC at their optimized heights increases dramatically with the FOV. Therefore, two LCWCs are enough for 1-arcmin FOV, but larger FOV will need more LCWCs.

Fig. 8. Influence of FOV on the optimal result of LCWC. (a) Reconstruction error JLC for two LCWCs as a function of FOV and conjugation height. Curves from the bottom to the top vary in θFOV from 0.2 arcmin to 2 arcmin with the increment of 0.2 arcmin. (b) Reconstruction error JLC at conjugation height of 10 km for two and three LCWCs.
3.4. Simulation of wavefront reconstruction for MCAO

Based on the results and discussion above, we use the optimal parameters as listed in Table 3 to have a performance testing. The performance is represented by Strehl ratio (SR). SR maps for the two MCAO systems are shown in Fig. 9.

Fig. 9. SR maps for two MCAO working with parameters listed in Table 3. (a) 1 arcmin FOV; (b) 2 arcmin FOV. The variances have been linearly interpolated on a regular grid of 317 evaluation directions. The green pentagrams sign the angle of the LGSs. The dotted circle is the range of FOV.
Table 3.

Parameters for MCAO.

.

The average SR will be 0.6 and 0.58 in the FOV. The two MCAOS both have good performance at most region of FOV. The SR varies in the observed target direction αt. In MCAO, we get higher SR of 0.74 at the direction of LGS. The minimum SR of 0.3 appears at the intersection of edge of FOV and the central line of two LGSs.

4. Conclusions

In this paper, a Zernike polynomial-based modal method was introduced to predict the higher-order error of MCAO across the field for arbitrary configuration of LGS and LCWC analytically. Then we obtained a phase reconstruction and compensation with high precisions by optimizing the location of LGSs and LCWCs conjugation heights. A demonstration of multi-LGS optimization process was given and discussed. We concluded that the optimal θg is linear to θfov and the optimal conjugated height of LCWCs is not affected by the FOV. Theoretical results indicated that our method would get the minimum wavefront reconstruction error for MCAO, and had higher imaging resolution at large FOVs. Additional experiments in our lab have been prepared to validate the method for the next-step research. Our method will benefit the design of MCAO in large aperture telescopes.

Reference
1Xuan LLi D YLiu Y G 2015 Chin. J. Liq. Cryst. Displays 30 (in Chinese)
2Booth M J 2014 Light-Sci. Appl. 3 e165
3Boyer CAdkins SAndersen D Ret al. 2014 Proc. SPIE 9148 91480X
4Bouchez A HActon D SBiasi RConan REspeland BEsposito SFilgueira J MGallieni DMcLeod B APinna ESantoro FTrancho GVan D M A 2014 Proc. SPIE 9148 91480W
5Ma J QLiu YChen J JLi B QChu J R 2014 Opt. Precis. Eng. 22 2047 (in Chinese)
6Liu R XZheng X LLi D YXia M LHu L FCao Z LMu Q QXuan L 2014 Chin. Phys. 23 94211
7Liu CM Q QHu L FCao Z LXuan L 2010 Chin. Phys. 19 64214
8Wang LLin X WHu WShao G HChen PLiang L JJin B BWu P HQian HLu Y NLiang XZheng Z GLu Y Q 2015 Light Sci. Appl. 4 e253
9Zhang Y CHu L FPeng Z HXuan L 2014 Chin. J. Liq. Cryst. Displays 29 709 (in Chinese)
10Chen Y CXuan LYao L SPeng Z H 2016 Chin. J. Liq. Cryst. Displays 31 62 (in Chinese)
11Yu YZhou XZhang Z D 2015 Chin. J. Liq. Cryst. Displays 30 213 (in Chinese)
12Wang Y CZhang Z DZhou X 2015 Chin. J. Liq. Cryst. Displays 30 949 (in Chinese)
13Min LYao L SWang Y KWang S XLiu Y GPeng Z HXuan L 2015 Chin. J. Liq. Cryst. Displays 30 32 (in Chinese)
14Wizinowich PActon DXa Set al. 2000 Publ. Astron. Soc. Pac. 112 315
15Zhang Z CYou ZChu D P 2014 Light Sci. Appl. 3 e213
16Johnston D CWelsh B M 1994 J. Opt. Soc. Am. 11 394
17Ellerbroek B L 1994 J. Opt. Soc. Am. 11 783
18Rigaut FNeichel BBoccas Md’Orgeville CVidal FVan D M AArriagada GFesquet VGalvez R LGausachs G 2014 Mon. Not. R. Astron. Soc. 437 2361
19Neichel BRigaut FVidal Fvan Dam M AGarrel VCarrasco E RPessev PWinge CBoccas Md’Orgeville C 2014 Mon. Not. R. Astron. Soc. 440 1002
20Gilles LEllerbroek B L 2008 J. Opt. Soc. Am. 25 2427
21Neichel BFusco TConan J M 2009 J. Opt. Soc. Am. 26 219
22Rosensteiner MRamlau R 2013 J. Opt. Soc. Am. 30 1680
23Fusco TConan J MMichau VMugnier L MRousset G 1999 Opt. Lett. 24 1472
24Ellerbroek B LGilles LVogel C R 2003 Appl. Opt. 42 4811
25Fusco TConan J MMichau VRousset GMugnier L M 2000 Proc. SPIE 4007 454
26Noll R J 1976 J. Opt. Soc. Am. 66 207
27Fan M WHuang L HLi MRao C H 2015 Opt. Precis. Eng. 23 2803 (in Chinese)
28Wang Y KHu L FWang C CWang S XXuan L 2015 Opt. Precis. Eng. 23 2203 (in Chinese)
29Tokovinin ALe Louarn MViard EHubin NConan R 2001 Astron. Astrophys. 378 710
30Ragazzoni RMarchetti ERigaut F1999Astron. Astrophys.342L53
31Whiteley M RRoggemann M CWelsh B M 1998 J. Opt. Soc. Am. 15 993
32Troxel S ERoggemann M CWelsh B M 1994 J. Opt. Soc. Am. 11 2100